library(dplyr)
library(ggplot2)
library(qqman)
library(ggrepel)
library(plotly)
library(manhattanly)

mhp_int <- function(data_loc) {
  dat <- read.table(data_loc, header = T) %>% 
    filter(!is.na(TWAS.P) & TWAS.P != 0) %>%
    rename(P = TWAS.P, SNP = BEST.GWAS.ID, GENE = ID) %>% 
    mutate(BP = (P0 + P1)/2, 
           LOGP = -log10(P)) %>% 
    select(SNP, GENE, CHR, BP, P, LOGP) %>% 
    group_by(CHR)
  manhattanly(dat, snp = "SNP", gene = "GENE",
              annotation1 = "P", annotation2 = "LOGP")
}
mhp_info <- function(data_loc, logp_bar = 3) {
  dat <- read.table(data_loc, header = T) %>% 
    filter(!is.na(TWAS.P) & TWAS.P != 0) %>%
    rename(P = TWAS.P, SNP = ID) %>% 
    mutate(BP = (P0 + P1)/2) %>% 
    select(SNP, CHR, BP, P) %>% 
    group_by(CHR)
  dat_plot <- dat %>% 
    summarise(chr_len = max(BP)) %>% 
    mutate(tot = cumsum(chr_len) - chr_len) %>% 
    select(-chr_len) %>% 
    left_join(dat, ., by=c("CHR"="CHR")) %>% 
    arrange(CHR, BP) %>%
    mutate(BPcum = BP + tot) %>%
    mutate(is_annotate = ifelse(-log10(P) > logp_bar, "yes", "no")) 
  axisdf <- dat_plot %>% 
    group_by(CHR) %>% 
    summarize(center = (max(BPcum) + min(BPcum)) / 2)
  ggplot(dat_plot, aes(x=BPcum, y=-log10(P))) +
    geom_point(aes(color=as.factor(CHR)), alpha=0.8, size=1.3) +
    scale_color_manual(values = rep(c("grey", "skyblue"), 22 )) +
    scale_x_continuous(label = axisdf$CHR, breaks= axisdf$center ) +
    scale_y_continuous(expand = c(0, 0) ) + 
    geom_label_repel(data=subset(dat_plot, is_annotate=="yes"), aes(label=SNP), size=2) +
    theme_bw() +
    theme( 
      legend.position="none",
      panel.border = element_blank(),
      panel.grid.major.x = element_blank(),
      panel.grid.minor.x = element_blank()
    )
}

ad_Brain_Amygdala

mhp_info("./Data/ad_Brain_Amygdala.dat")

mhp_int("./Data/ad_Brain_Amygdala.dat")

ad_Brain_Cerebellum

mhp_info("./Data/ad_Brain_Cerebellum.dat")

mhp_int("./Data/ad_Brain_Cerebellum.dat")

ad_Brain_Cortex

mhp_info("./Data/ad_Brain_Cortex.dat")

mhp_int("./Data/ad_Brain_Cortex.dat")

ad_Brain_Hippocampus

mhp_info("./Data/ad_Brain_Hippocampus.dat")

mhp_int("./Data/ad_Brain_Hippocampus.dat")

t2d_Adipose_Subcutaneous

mhp_info("./Data/t2d_Adipose_Subcutaneous.dat")

mhp_int("./Data/t2d_Adipose_Subcutaneous.dat")

t2d_Adipose_Visceral_Omentum

mhp_info("./Data/t2d_Adipose_Visceral_Omentum.dat")

mhp_int("./Data/t2d_Adipose_Visceral_Omentum.dat")

t2d_Pancreas

mhp_info("./Data/t2d_Pancreas.dat")

mhp_int("./Data/t2d_Pancreas.dat")

t2d_Whole_Blood

mhp_info("./Data/t2d_Whole_Blood.dat")

mhp_int("./Data/t2d_Whole_Blood.dat")